INFLUENCE OF REFRACTORY PERIODS IN THE HOPFIELD MODEL 
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I. INTRODUCTION 



During the last years some biological features of real neurons have been incorporated into the Hopfield model 
|E| in order to make it more realistic and trying to improve its performance. Suitable modifications of the original 
model taking into account biological ingredients such as thermal noise, dilution, asymmetry, dynamical delays, among 
others, have been vastly analized in the literature ]2]-p^]. Although they usually deteriorate the retrieval ability, it 
has been shown they enable the implementation of new tasks, such as recognition of temporal sequences [fz| — |l~2|| and 
categorization jT3-[l3]. 

One crucial biological element originally absence in the Hopfield model is the so called refractory period [[L6|- 
Real neurons take about 1-2 milliseconds to complete a cycle from the emission of a spike in the pre-sinaptic neuron 
to the emission of a spike in the post-sinaptic neuron. After this time, the neuron need again about 2 milliseconds to 
recover, and during this time, called absolute refractory period (ARP), it is insensitive to afferent activiy (i.e., it cannot 
emit a second spike, no matter how large the post-synaptic potential (PSP) may be). Following this short ARP, the 
neuron enters in a new regime of about 5-7 milliseconds, in which it partially recovers the capacity of emitting spikes, 
but now with a greater excitation threshold which decreases with time. This is called the relative refractory period 
(RRP). Following this somewhat longer RRP, the threshold tends to return to its rest value and the neuron can fire 
again with typical intra- network potentials. 

The simplest way one can introduce these periods into the dynamics of the Hopfield model is by means of a time 
dependent threshold acting as an external field, which depends on the recent history of the neuron. Since we want 
this threshold to mimic the effects of fatigue observed in real neurons , it should act only after the cell has emitted 
an electric signal. So, we expect that the threshold depends on the mean activity of the neuron in the previous time. 
The main effect on the dynamics of the model is to introduce a tendency to destabilize the fixed point attractors, 
allowing the appearance of oscillatory behaviors. In the last years different threshold functions have been studied 
p| , ^9[ -21|, showing that they enable the system to wander through the phase space, eventually visiting different basins 
of attraction and simulating the process by which the brain recognizes temporal sequences of patterns . On the other 
hand, oscillating and chaotic trajectories in the phase space seem to be more realistic than fixed points attractors 
from a biological point of view (see Jl^] and references therein). 

In this work we analyze, using a mean field approach and through numerical simulations, the behavior of the 
Hopfield model for associative memory when the effect of these refractory periods are taken into account in the 
dynamics of the system. Instead of considering a fatigue like threshold function that would depend on the large term 
history of the neuron [jll 20 1, we introduce a threshold that depends only on the state of the neuron in the previous 
time, i.e., it is activated only when the neuron fires a spike. In the section |H], we introduce the model and describe 
how the refractory periods are incorporated into its dynamics. In section M , we obtain an equation for the value 
of the superposition between the state of the system and one of the memories (which is only valid for fixed points 



dynamics), from which we can study the retrieval properties of the model in this region. In section IV, we obtain 



a complete phase diagram and identified the regions of fixed points, cyclic orbits and chaotic orbits. We have used 
a synchronous parallel updating, which allows an efficient use of modern parallel-processing computers. Finally, in 
section [v], we discuss the main results. 



II. THE MODEL 

As in the Little p2| and Hopfield models we consider a network of N binary neurons, each one modeled by an Ising 
variable Si which take the values {—1, +1}, representing the passive and active states, respectively. In order to take 
into account the effect of the refractory period in the neuron i we add a threshold that depends on the time, but only 
through the value of the state Si(t) of the neuron i. So the post-synaptic potential at time t is given by: 

hi(t) = hf(Jt)-j{l + Si{t)) , (1) 
where is the usual Hopfield post-synaptic potential: 

N 

hf(t)=J2 J t3 S,(t) . (2) 

Here J, t j is the Hopfield synaptic matrix connecting the pre and post-synaptic neurons j and i and whose elements 
have the form: 
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jU=l 

The £^ are random independent variables which take the values ±1 with the same probability and the N-bits words 
{(j, £2 > • ■ • j Cjvi stand for the p stored configurations (fJ, = 1,2,... ,p). The dynamics of the network is governed by 
a Monte Carlo heat bath dynamics: 

-1 with probability ( 1 + exp(— ) ) 

-1 with probability ^1 + exp(+ ) j 

where all the neurons are updated simultaneously (like in the Little model). The parameter T measures the noise 
level of the net and in the noiseless limit (T = 0) we recover the deterministic dynamics: 

Si(t + 1) = Sign (hf(t) - ^(1 + (5) 

From this expression we can easily understand the effect of this extra field: if the neuron i fires a spike at time t 
(Si(t) = +1), it will requires an extra contribution A to the PSP in order to fire again. On the other hand, if this 
neuron was at rest at time t (Si(t) = —1), then it will work like an usual Hopfield neuron. Observe that this model 
does not distinguish between absolute and relative periods neither includes any fatigue like effect (long time history) . 
As usual, we will characterize the recognition ability by calculating the long time behavior of the overlap between 
the state of the system {Si(t)} and the stored patterns, defined as: 



1 N 

"VW"E <S(«)>r (6) 



where (. . .) means a thermal average at temperature T. We say that the system recognizes a pattern every time it 
evolves to an attractor for which only one overlap is non-zero and all the others vanish as (0(1/ y/(N))). The two 
relevant parameters in our model are then A and a (the ratio between the number of stored patterns (p) and the 
total number of neurons of the network (N)). In the following sections we analyze the behavior of the model on the 
(A, a) plane. 



III. FIXED POINT EQUATION 

Following the statistical method developed by Geszti |^3| (see also p^|27|), we give in this section a heuristic 
derivation of the critical capacity a c as a function of the parameter A for the stochastic version of the model. By 
taking the limit T->0we obtain a noiseless phase diagram in the (A, a) plane which will be compared with numerical 
simulations in the next section. 

Let us suppose that the initial state of the system is such that mi = m is the only macroscopically non-zero overlap 
and so = 0(1/ y/(N)) for any [ij^l. Furthermore, we will assume that although the threshold tends to destabilize 
the fixed point attractors, its effect is not strong enough to anable the system to visit different basins of attractors. 
So, since initially only the first overlap was non zero, let us suppose that this will be valid for any time t. This a priori 
assumption will be justified in the next section by the numerical simulation, where we will find that in the region 
where the system recognizes (that is, where m = mi(t — ► 00) ^ 0) the dynamics of the model is dominated by fixed 
point attractors. 

We then start considering the overlap between the state of the system and the first pattern, that can be rewritten 
as: 

Since we are storing an extensive number of pattern, we cannot neglect any more the effect of the others (p — 1) 
overlaps: 



3 



m + ere 1 ™, + E ti$ m n - # t C 1 + s *) ] ] ( 8 ) 



In order to make an self-consistent treatment for the overlap m we need to introduce two other parameters, namely: 

jv 

i 

1 



= ^E<^> 2 w 

i=l 
1 P 

r = ~Y. m l ( 10 ) 



where q is the Edwards- Anderson order parameter and r is indentified as the mean square overlap of the system 
configuration with the nonretrieved patterns p5| . 

After some standard calculations we get the following set of equation for the values of to, q and r m i/ie attractor. 

m =\ J Dz ( tanh (P L f) + tanh (P L 7)) ( n ) 
o=\ J Dz {tanh 2 (f3L+) + tanh 2 (/?£")) (12) 



'1 



(13) 



where Lf = (1 — y)m ± y + yfarz and 

(iz exp —z 2 j1 



Dz 



'2tt 



Notice that for the particular case A = we recover the equations obtained for the Hopfield model pq] which also 
agree with those obtained by Amit et al 0] through a thermodynamical mean-field study (which unlike this method 
requires the use of the replica trick). 

We start analyzing the noiseless case (T = 0) for which we have performed numerical simulations. In this limit our 
equations take the following form: 

.A\ mJ .i\ -i /f1_A^_A\ 

(14) 





(l-C) 



(16) 



In Fig. 1 we display to as function of a for several values of A. For any value of A < A c = 1 there always exists a 
critical value a c below which the system recovers the stored patterns with a non-zero fraction of errors e. At a c (A) the 
systems undergoes a discontinues transition from the retrieval phase (in which the dynamics is governed by the fixed 
point attractors) to a non-retrieval phase where our analytical approach is no longer valid, since the self-consistent 
equation does not predict a fixed point attractor (which was our original assumption). Observe that a c decreases as 
A increases. As a — > the fraction of errors at the transition e c = i(l — to) goes to accordingly to the following 
expression: 




1 exp{- 



(l-A) 2 



e ^ )+ 1-A ) (17) 

We have also analyzed the fixed point equations in the presence of noise. In Fig. 2 we present the (T, a) phase 
diagram for different values of A. For A = we recover the phase diagram obtained in pj. Along the lines T c (a) 
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the system undergoes a discontinues transition from the retrieval phase (below) to the non-retrieval phase (above). 
Notice that the recognition phase decreases as A increases; i.e., the main effect of introducing this refractory periods 
seems to be a degradation of the retrieval properties of the model. In Fig. 3 we present the critical line T versus 
A for a — 0. For A < 0.611 the system undergoes a second order transition while for A > 0.611 the transition is 
discontinues (the point (T, A) ~ (0.46,0.611) separates both lines). In the inset we show the behavior the retrieval 
overlap around the critical point as function of T. 

IV. NUMERICAL SIMULATION 

In this section we present a numerical study of both recognition ability and dynamical properties of the model at 
T = and compare it with the analitical results obtained in the previous section. The simulations were performed 
on systems of N = 800, 1600 and 3200 neurons and the network was updated synchronously. Setting the initial 
configuration as the first stored pattern, we let the system evolve until it reaches the attractor. 

In order to characterize the dynamical behavior we first determined whether the system was in a periodic orbit or 
not, by waiting until it returned to a given configuration that was stored after a transient. Depending on the value 
of the parameters and on the size of the system it could also happen that the system did not return to the initial 
configuration after a given period of time (typically 100 Monte Carlo Steps). In such cases, we said that the system 
follows a chaotic orbit, although we have not performed a through analysis in order to determine whether these were 
really chaotic orbits or orbits with large periods. 

To analyze the recognition ability we calculated for each sample a temporal average between the stored patterns 
and the state of the system in the attractor. If the system reached a cyclic orbit of period t c , we measured (in the 
attractor) the following quantity: 

^ to+t c 

"V = ^ E m * (*) ( 18 ) 

c t=to 

Since the initial state was chosen to be always the first memory, we say that the network recognizes when 

m = mi ~ 1 ~ 0(l/y/(N)), for n > 1 (19) 

In order to make a configurational average of to, for any value of the parameters we repeated this procedure over 100 
different samples using different memories, initial configurations and random number sequences. To characterize the 
dynamical behavior we present the frequency with which each kind of attractor appears and also the mean activity, 
defined as the average number of active neurons, in the attractor. 

1 N 

=2ArE< 1 + fl W (20) 

1=1 

In Fig. 4 we display the phase diagram A vs. a for N = 3200. For A = the system presents only fixed points 
(FP). For fixed a, as A increases we found that: 

1. for low values of A the dynamics is governed only by fixed points attractors. The full circle indicates where this 
kind of behavior disappears; 

2. the region between the two full triangles indicate the region where cycles of order two (C2) appear; 

3. the hollow circle indicates the value of A above which chaotic orbit (Ch) emerges. 

Observe that there are many region of coexistence of attractors. In fact, between the C2 and the Ch we have also 
found cyclic orbits (OC) of order greater than two, but they are not indicated in the diagram. 

Independently of the dynamical behavior, we have also studied the critical recognition capacity. The dashed line 
separates the recognition phase (below) from the non-retrieval phase (above) obtained numerically and the full line 
corresponds to the analytical results obtained in the previous section. The simulation curve fits very well the analytical 
result only for small values of a. 

In order to understand why the analytical and the numerical curves do not agree, we have carefully analyzed the 
behavior of the system along two cuts with fixed a, namely 0.01 and 0.04. In Fig. 5 we plot both to (top) and the 
frequency with which each kind of orbits appears (bottom) as a function of A. The first thing we note is that the 
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FP region coincides with the retrieval phase, and that the C2 region corresponds to the non-retrieval phase. In such 
cases, where the systems only recognizes through FP, the analytical curve predicts very well the transition. On the 
other hand, in Fig. 6 we present the same curves for a = 0.04. Notice now that the recognition phase presents two 
different dynamical behaviors: for small values of A the system evolves to FP while for intermediate values it goes to 
C2. Unlike the a = 0.01 case, now the theoretical curve does not predict correctly the retrieval non-retrieval phase 
transition, but the FP to C2 transition. In both cases we have studied the finite-size effects by working with three 
different sizes, namely, N = 800, 1600 and 3200. In Figs. 5 and 6 we present the overlaps as function of A for all 
these system sizes. Note that as N increases the numerical simulation tends to display a more abrupt decay of m at 
the transition, resembling the first order transition found in the analytical calculation. 

Finally, in Fig. 7 we show the mean activity as function of A for a = 0.01 and 0.04. We can notice that where 
there are fixed points and periodic orbits with recognition within, the mean activity remains around the value ~ 0.5 
(random variables) and it only decreases in the transition to non-retrieval phase. This shows that the parameter A 
not only damages the recognition ability but also destabilizes the tendency of the system to evolve to fixed point 
attractors, allowing the appearence of more complicated retrieval attractors. 

V. CONCLUSIONS 

In this work we study analytically and through of numerical simulations a model for associative memory where we 
have incorporated in the dynamics of the network a new kind of threshold that simulate the effect of the refractory 
period. The main result is that the parameter A that activates this threshold yields to the appearing of Chaotic and 
periodic attractors. Nevertheless, the system seems to recognizes only through fixed point and cycles of order two. 
Only in a small region the system recognizes with higher order cycles and with chaotic trajectories, but this behavior 
appears just in the boundary between the retrieval and the non-retrieval phases. It would be interesting to make a 
more detailed study to elucidate whether this kind of trajectories are due to finite size effects or not. As much as we 
could see, as N increases they do not seem to dissapear, so we suspect that they will exist also in the thermodynamical 
limit. 

In the recognition phase (small values of A), the PSP is strong enough to drive the system to stable attractors, FP 
and periodic orbits, where the average overlap in each regime is of the order 1. For large values of A the performance 
is drastically damaged, and in these regions the dynamics is dominated by very large cycles or chaotic trajectories. 

The numerical simulation fits very well the analytical results only for small values of a, where the transition occur 
from fixed point FP to cycle order two C2. Actually, the analytical curve seems to fit only the line where the fixed 
point behavior disappears. We also observe that in the transition the mean activity decreases with the increase of A. 
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CAPTIONS FOR FIGURES 



Figure 1. Plot of m versus a at T = for different values of A. At a c (A) the system undergoes a discontinuous 
transition from the recognition phase to non-retrieval phase. 

Figure 2. Phase diagram T versus a for A = 0, 0.2, 0.4 and 0.6. Below of the critical lines the system recognizes 
with fixed point and the transition to non-retrieval phase is discontinues. 

Figure 3. The critical line T = /(A) for a = 0. For A < 0.611 the transition is of second order (full line) while 
for A > 0.611 the transition is discontinues (dashed line). (T, A) ~ (0.46,0.611) is a critical point. 

Figure 4. The numerical phase diagram A versus a at T — and N = 3200, showing the regions FP (below full 
circles), periodic (between the two full triangles) and Ch (above hollow circles). The simulation (dashed line) and 
analytical (full line) curves separetes the recognition phase (below) from non- retrieval phase (above). 

Figure 5. Plot of m (top) and of the frequency (bottom) in the which each kind of orbits appears as a function of 
A for a — 0.01. The full line corresponds to the analytical curve. 

Figure 6. Plot of m (top) and of the frequency (bottom) in the which each kind of orbits appears as a function of 
A for a — 0.04. The full line corresponds to the analytical curve. 

Figure 7. The mean activity vs. A for a = 0.01 and a = 0.04. 
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